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Motivated by the recent demonstration of its use as a tool for the detection and characterization of 
phase-shape correlations in multivariate time series, we show that eigenvalue decomposition can also 
be applied to a matrix of indices of bivariate phase synchronization strength. The resulting method 
is able to identify clusters of synchronized oscillators, and to quantify their strength as well as the 
degree of involvement of an oscillator in a cluster. Since for the case of a single cluster the method 
gives similar results as our previous approach, it can be seen as a generalized Synchronization Cluster 
Analysis, extending its field of application to more complex situations. The performance of the method 
is tested by applying it to simulation data. 
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1. Introduction 

The eigenvalue decomposition of the equal-time correlation matrix of a set of signals is 
one of the standard tools of multivariate data analysis (cf. lAndersonl (2003 1). Recently, 
Miille ref al. 1 1|2005 1 demonstrated the usefulness of the eigenvalue decomposition of the 
correlation matrix specifically as a tool for the detection of phase-shape correlations in 
multivariate data sets. They showed that changes in the degree of synchronization in all or a 
subset of signals are reflected in coordinated changes in the highest and lowest eigenvalues, 
and that information about the channels involved and the type of their interaction can be 
obtained from the corresponding eigenvectors. 

While the correlation of time series may indicate the synchronization of the oscillators 
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they are obtained from, the physical concept of synchronization refers specifically to the 
adjustment of the rhythms of oscillators, i.e. to the relative dynamics of their phases rather 
than their amplitudes [Pikovsky et al. (2001) |. Moreover, there is a regime in the dynamics 



of coupled chaotic oscillators in which the phase difference is bounded while the ampli- 
tudes remain uncorrelated [Rosenblum et al. (1996)], called phase synchronization. In this 
paper, we show that in order to focus the analysis on synchronization relations, it is possible 
to replace the matrix of correlation coefficients with a matrix of indices of bivariate phase 
synchronization strength. Combined with an additional step of sorting signals into groups, 
eigenvalue decomposition can operate as a Synchronization Cluster Analysis, generalizing 
the previous approach of Allefeld and Kurths [2004|. 

2. Eigenvalue Decomposition of the Synchronization Matrix 

The correlation matrix C of a set of data channels x,, i = 1 . . . N, consists of the correlation 
coefficients C (/ 6 [-1; 1] between channels. Its eigenvalues A^ and eigenvectors % are 
defined by the equation 

C v k = A k v k , (1) 

which in general has N different solutions, k = 1 . . . N. In the following we assume that the 
eigenvectors are normalized, \v k \ = 1, and the solutions have been sorted according to the 
eigenvalues, A\ < A% < . . . < An- 

The eigenvectors and -values of C are real-valued, and the eigenvalues are non-negative. 
Because a matrix becomes the diagonal matrix of its eigenvalues by being transformed into 
the basis of its eigenvectors and the trace of a matrix is invariant under such a transform, for 
every correlation matrix holds 2 A k — tr(C) = N. In the uncorrelated case this is trivially 
fulfilled by A k = 1 for all k. With a deviation from this, each increase of an eigenvalue 
above 1 has to be compensated for by at least one other eigenvalue becoming smaller than 
1, such that this value gives a natural distinction between "large" and "small" eigenvalues. 

The quantification of phase synchronization is based on the instantaneous phase (pi of 
each oscillator i = 1 . . . N. How these phases are determined in the special case is not 
important here; if the given data are time series, the standard approach is the Hilbert trans- 
form for narrowband data, or the Morlet wavelet transform for broadband signals (for a 
discussion, see|Pikovsky et a/. | II 19971 or Allefeld [2004|, Sec. 3.2). The statistical strength 



of phase synchronization of two oscillators i and j can then be defined as the "peakedness" 
of the distribution of the phase difference <pj - <pi; here we use the measure 



(2) 



1 

-/.exp(i(0j/-0a)) , 
n 1~t V ' 

where 1=1 ...n enumerates the realizations in the given sample. For the continuum from 
no to perfect phase synchronization, this measure takes on values from to 1. R can be 
seen as the modulus of the complex correlation coefficient of signals x-, = exp(i </>,), and its 
decomposition shares the properties given above for C. 

Since R is a nonlinear measure, in this case the eigenvalue decomposition can no longer 
be interpreted with regard to a linear transform of data channels into source channels. But 
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Fig. 1. Left: Synchronization matrix consisting of three clusters of oscillators. Right: Its eigenvectors and -values. 
Eigenvectors corresponding to eigenvalues > 1 describe the cluster structure. 
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Fig. 2. Synchronization matrix consisting of three clusters with additional inter-cluster synchronization. The clus- 
ter structure appears in the eigenvectors for > 1, but only in different superpositions. 



still the result of the decomposition can be used as a means to analyze the structure of 
synchronization relations. We will demonstrate this with two basic, artificially constructed 
examples of synchronization matrices. 

Figure Q] shows the synchronization matrix for a system consisting of three clusters 
of synchronized oscillators (with no synchronization between clusters and a different de- 
gree of involvement of each oscillator in its cluster) along with the result of the eigenvalue 
decomposition. There are three eigenvalues larger than one, and the corresponding eigen- 
vectors describe clearly the extent of the clusters as well as the degree of involvement of 
the individual oscillators. There is also a correspondence between the size and strength of 
internal synchronization of the clusters and the three eigenvalues. In contrast, the remaining 
eigenvectors seem not to contribute to the description of the synchronization clusters. The 
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original participation indices trimmed matrix clusters 
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Fig. 3. Synchronization Cluster Analysis. Left: Participation indices corresponding to the > 1 for the synchro- 
nization matrix of Fig. [2] Oscillators are attributed to that cluster for which its participation index is maximal: 
blue, green, or red. Center: Synchronization matrix with removed inter-cluster synchronization. Right: Result of 
the eigenvalue decomposition of the trimmed matrix; shown are the participation indices and cluster strengths. 

interpretation of the eigenvalue decomposition of R therefore has these aspects: 1 ) Synchro- 
nization clusters are identified by eigenvalues Ak> 1 . The eigenvalues themselves quantify 
the strength of the clusters. 2) For each cluster, the corresponding eigenvector describes 
its internal structure. Because v\ = 1, the index v\ quantifies the relative involvement 
of channel (oscillator) i in cluster k. 3) Combining both, the "absolute" involvement of 
channel i in cluster k can be quantified by the participation index /Uvl. 

Figure [2] gives the result for a matrix consisting of the same three clusters, but with 
additional inter-cluster synchronization. Because of the coupling between them, the three 
clusters no longer appear in separate components of its eigenvalue decomposition. There 
are still three X% > 1, but the eigenvectors consist of superpositions of the clusters. To ac- 
count for this, the oscillators belonging to the three clusters have to be identified explicitely. 

This can be done by means of the participation indices; the procedure is illustrated in 
Fig. [3] Each oscillator is attributed to that cluster for which its participation is maximal. 
In this way, the three clusters consisting of oscillators #1-5 (blue), 6-9 (green), and 10-12 
(red) are correctly identified. In a second step, all of the indices for inter-cluster synchro- 
nizations are set to zero, and the eigenvalue decomposition is repeated on the trimmed 
matrix. For the result of this decomposition, the interpretation given above is valid again. 

3. Generalized Synchronization Cluster Analysis 

The procedure described in the last section is an approach to synchronization cluster anal- 
ysis. The method identifies clusters of synchronized oscillators and quantifies the strength 
of the clusters as well as the degree of involvement of each oscillator in its cluster. In a 
previous paper, Alle feld and Kurtrisl 020041 introduced another approach that was limited 
to a single cluster of synchronization; it assumed that all of the oscillators belong to the 
same cluster, and focused on the quantification of the degree of oscillator participation. 
The algorithm derived from the observation that under relatively general conditions the 
synchronization indices within a cluster can be written as the product of factors R,c, 



Ry = RcRjc for i * j (R« = 1), 



(3) 
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Fig. 4. Relation between the result of the single-cluster analysis (horizontal scale: estimate of R^ c ) and the parti- 
cipation indices given by eigenvalue decomposition (vertical scale: Afjv 2 N ) in 20 simulation runs. 



which can be interpreted as the synchronization strength between oscillator i and the cluster 
itself, its "to-cluster synchronization strength". 

We investigated the relationship between the results of the two methods by applying 
them to simulation data that conforms to the presupposition of a single factorizable syn- 
chronization cluster. For N = 20 oscillators, the R,c were drawn randomly from the uniform 
distribution over [0; 1]; then the R/y were calculated based on 100 samples of the phase 
difference from a wrapped normal distribution (cf. Allefeld andKurths [ 2004-1). Figure [4] 
shows the result for the relation between the output of the earlier method and the participa- 
tion indices for the strongest cluster. Though the mathematical background of both methods 
is clearly different, the plot shows an almost functional dependency, which can be roughly 
described by R^-. = /tjvvjL. The region of zero participation indices for R^-. below about 0.2 
comes from the attribution of weakly synchronized oscillators to another cluster by the new 
method; in this range the observed value can no longer be reliably distinguished from no 
synchronization for the given sample size. Since there is a clear relationship between the 
results of both methods in the case where the assumptions of the earlier method hold, the 
cluster analysis based on the eigenvalue decomposition of the synchronization matrix can 
be seen as a generalization of our previous approach to synchronization cluster analysis. 



4. Application to Simulated Phase Synchronization 

To check the performance of the new method, we applied it to data obtained from the nu- 
merical simulation of a system that is known to exhibit clusters of phase synchronization, 
previously investigated by Osipov et al. 1 1997 1 (see Eq. 1 & Sec. IV A). The system con- 
sists of a chain of 20 phase-coherent Rossler oscillators with a diffusive coupling in the 
y-component of strength e = 0.007. The natural frequencies of the oscillators increase lin- 
early along the chain, in the range a> — 1 ... 1 .004. Phases were defined as the angle in 
the x/y-plane. To also test for the influence of small sample size and noise, we addition- 
ally applied the algorithm to the synchronization matrix obtained from data reduced to 100 
approximately independent samples with 50% white noise added to the time series data. 

The result is shown in Fig. [5] The analysis identifies two clusters of synchronized os- 
cillators (#1-16 & 17-20). There is an almost constant high participation of oscillators 
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Fig. 5. Analysis results (participation indices and cluster strengths) for a chain of nonidentical Rossler oscillators. 
Left: Based on original simulation data. Right: Data reduced to 100 samples, plus 50% measurement noise. 

in their cluster, resulting in cluster strengths close to the number of involved oscillators. 
Decreased participation occurs near the border between the clusters, where the coupling 
along the chain pulls oscillators away from the common dynamics of their group. The re- 
sult for reduced sample size and measurement noise is very similar to the original. Though 
participation indices and cluster strengths are slightly smaller due to the noise, the same 
two clusters are clearly identified. This result indicates that the cluster analysis based on 
eigenvalue decomposition is relatively robust against small sample size and noise. 

Conclusion and Acknowledgements 

We introduced a new approach to synchronization cluster analysis based on eigenvalue 
decomposition. The new method can be seen as a generalization of our previous approach, 
extending its field of application to situations including multiple clusters. The algorithm 
was tested on simulation data and shown to be robust against small sample size and noise. 
Future work will be on the use of the method to investigate synchronization patterns in EEG 
data. A Matlab implementation of the algorithm can be obtained from the corresponding 
author. — This work has been supported by Deutsche Forschungsgemeinschaft. 
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